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■ 

■ Using techniques from singular perturbation theory, we explicitly calculate the cosmological evo- 
' lution in a class of modified gravity models. By considering the (m)CDTT model, which aims to 

, ^ , explain the current acceleration of the universe with a modification of gravity, we show that Einstein 

■ evolution can be recovered for most of cosmic history in at least one f{R) model. We show that a 
' standard epoch of matter domination can be obtained in the mCDTT model, providing a sufficiently 
I long epoch to satisfy observations. We note that the additional inverse term will not significantly al- 
. ter standard evolution until today and that the solution lies well within present constraints from Big 
• Bang Nucleosynthesis. For the CDTT model, we analyse the "recent radiation epoch" behaviour 
; (a (X t^/^) found by previous authors. We finally generalise our findings to the class of inverse 

1 I , power-law models. Even in this class of models, we expect a standard cosmological evolution, with 
i-^ ' a sufficient matter domination era, although the sign of the additional term is crucial. 

6 ■ 

^ ; I. INTRODUCTION 

, The observation that the universe appears to be accelerating at present times has caused one of the greatest problems 
^ to modern cosmology. High precision data from Type la supernovae, the Cosmic Microwave Background and large 
scale structure seem to hint that the universe is presently dominated by an unknown form of energy, dubbed dark 
energy [1, 2, 3, 4, 5, 6]. The problem resides in the lack of good physical motivation for this energy density. One 
' obvious contender for the role of dark energy is Einstein's cosmological constant, but particle physics fails to predict 
. the correct density. Other candidates include scalar fields or extra dimensions. We refer the reader to [7] for a review. 
^-H ■ Recently, a modification of General Relativity itself was suggested to explain this accelerating universe [8, 9, 10, 
11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21]. For reviews see e.g. [22, 23, 24]. The postulation is that in which 
the Ricci scalar of the Einstein-Hilbert lagrangian is replaced by an arbitrary function of the Ricci scalar, f{R). 
The most simple example includes inverse powers of R: f{R) ~ R — a/R^ [15]. While this modified action was 
shown to give an accelerated attractor, due to the equivalence to scalar-tensor theories with a vanishing Brans-Dicke 
parameter[25, 26, 27], solar system experiments rule these models out {ujbd > 40,000 [28]). Attempts have been 
made to overcome this problem [17], but solar system gravity still raises doubt as to these models [29]. 
' These models have raised much recent interest, due to their perhaps simple nature. The presence of ghosts and 
stabilities have been studied [30, 31, 32], as well as the implications to the CMB, large scale structure, solar system 
gravity and gravitational waves [33, 34, 35, 36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49]. 

Recently, it has been claimed that, in all theories which behave as a power of R at large or small R, standard 
cosmological evolution can not be obtained [50, 51]. Further, it is claimed that no f{R) will lead to an acceptable 
matter era [52]. Specifically, during the matter phase, the authors find that in these theories the scale factor evolves 
as a oc t^/^ and not as t^^^ as required. In order for a modified theory of gravity to be convincing, we require a 
close-to-Einstein evolution in the past, with a standard, sufficiently-long period of matter domination. We note that 
several models have been proposed which lead to matter dominance before acceleration [41, 53, 54] and that this issue 
is a matter of debate in the literature. 
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If the modification to gravity exists as a perturbative term in addition to the Einstein-Hilbert terms, then it is the 
beUef of the authors of this paper that Einstein cosmology should be recovered in the limit that the perturbation 
term disappears. This system has great similarity with perturbation theory [55]. This paper studies the asymptotic 
expansion of the (m)CDTT model [15], to show that standard cosmological evolution can be produced in at least one 
model. Specifically, we aim to recover an acceptable matter era and show that there is at least a class of f{R) for 
which this is possible. For the purposes of this paper, we do not concern ourselves with the constraints from solar 
system experiment. 

The paper is laid out as follows. In Section II, we define the theoretical framework and give some useful formulae. 
We state the simplest example in Section III and derive the asymptotic solution in Section IV. The relevance of this 
solution to cosmology is highlighted in Section V, where we comment on the era of matter domination, as well as Big 
Bang Nucleosynthesis. After generalising our specific example to the class of inverse power-law models in Section VI, 
we conclude in Section VII. We discuss the eigenmode analysis in the Appendix. 



II. THEORETICAL FRAMEWORK 



For convenience, we consider a class of modified gravity in which we add a perturbative function, f{R), to the 
Einstein-Hilbert action. We consider the action 



Jd^xV^g(^^±l^+C^y (1) 



where R is the Ricci scalar and Cm is the matter Lagrangian. The field equations, using the metric approach, can be 
derived from this action, 

R,,. - ^Qi^uR + JrR^c. - ({ - 9^.u - ^i^^.Ir = k'T^. (2) 

where i?^,^ and are the Ricci and stress-energy tensors respectively and we define fn = df /dR. Later, we will also 
require fuR = d^f/dR"^. Although this theory can be rewritten as a scalar tensor theory (e.g. [17, 56]), we shall not 
use the conformal transformation; we complete all our calculations in the Jordan frame, given by the action above. 

In a flat Friedmann-Robcrtson- Walker (FRW) background, the Ricci scalar can be written in terms of the scale 
factor, a{t) and the Hubble parameter, H = a/a. An overdot denotes differentiation with respect to cosmic time. 
Specifically, 



i2 

^a d 

In this modified theory, the Friedmann equation becomes 



3i?' ~ 'dfR [h + H^)+^ + ISH/rr i^H + AHH) = n'p. (4) 

As expected, when f{R) = 0, we recover Einstein cosmology, = K^p. The continuity equation remains the same 
as in Einstein gravity: 

PT + 3H {pT + pt) = 0, (5) 

where pr and pr are the total energy density and pressure respectively. These can be separated into components of 
radiation and matter, for example pt = Pr + Pm- 

We assume the parameter change as in CDDETT [18]: 

X = -H{t) y = H{t) 

so that 

H = -y^ = -— + —— ^^ Qx^{2v-1) 

dx dx ' V 

As in CDDETT, we also define a function, v{x), which can be related to the power of a power-law solution, a{t) oc t^: 

v{x) = — - 



x^ 



y 
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Then, power-law solutions can be identified simply, since v{x) p = constant. In standard Einstein cosmology, 
solutions for v{x) are expected to plateau at w = i and u = |, which correspond to epochs of radiation and matter 
domination respectively. We also note that an accelerating phase must have v{x) > 1. Using these new variables, it 
is then possible to rewrite Eqn. (4) in terms of x,v and dv/dx. This equation is often simpler to solve, as found in 
CDDETT. 

One most natural way of maintaining Einstein cosmology in the past whilst also explaining late-time acceleration 
is to consider a perturbative function, f{R), (perhaps an inverse power- law) which incorporates some perturbation 
parameter. It will be convenient later to separate this small parameter from the remaining function: 

/(i?) = ef{R). (7) 

In terms of this new function, f{R), (and using the new variables, x and v) the Friedmann equation takes the form: 

3x^ — k^Pt = es(x, v)^ + et(x, v) (8) 
dx 



where 



3x^{v-l)~ f 36a;^(2t;- 1) - 

/H — 77 H n J. 



S = o-fRR 

yd 

Eqn. (8) can be rearranged to give the useful formula: 



I8x^ ~„ 



RR 



^dv_ ^ ^ t{x,v) ^ [3x^ - k'^pt] 
dx s{x,v) s{x,v) 



III. EXAMPLE: THE CDTT MODEL 



One of the simplest modifications which aims to give late-time acceleration was proposed in [15] (the CDTT model): 

f(^) = ^^- (10) 

In the original paper, cr = — 1. Wc will also consider the mCDTT model, where a = I. We do not wish to advocate 
this model as a realistic modification of gravity, due to the impossibility to reconcile this action with solar system 
gravity [25]. Additionally, in order to justify the observation of a presently accelerating universe (the coincidence 
problem), we require p, Hq, which is fine-tuned. Ignoring these problems, we use this simple model to demonstrate 
the techniques. 

From Eqn. (4), the dynamical equations in such an action are 

4 

3H^ + ( 2HH + IbH^H + 2H^ + 6H^) = K^p. 

i2(if + 2/72)3 V ; ^ 

We always consider the general case with both radiation and matter, such that p = pr + Pm- In the terminology of 
Eqn. (7), e = /i* and f{R) = al/R. We therefore obtain the system: 

^4 dt; ^ ^4 3(l-2^)(2-.).; _ ^^^^^^ _ ^ _ ^^^^ 



dx 2x 

dpr Avpr 

dx X 

dpm ^VPn 



dx 

Trivially, when /U = 0, we have: 



(12) 
(13) 



K'^{Pr+ pm)-^x'^ = or {2v - 1) = 



We ignore the solution, v = 1/2, which is picked up due to rearrangement of the Friedmann equation, Eqn. (11). 
With only radiation present {pm = 0) and using Eqn. (12), v = 1/2. With only matter (p,. = 0), v = 2/3. We, of 
course, wish to study the asymptotical analysis, when p j^O, which will be the focus of the next section. 
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IV. ASYMPTOTIC ANALYSIS OF THE CDTT MODEL 



We now study the CDTT model in the case when 7^ 0, but is a perturbative parameter. We wish to find the 
asymptotic sohition as both /i ^ (to compare to Einstein's solution in the past) and x = 0{^i) (with the aim to find 
late-time acceleration). In the limit /U ^ 0, (11)-(13) gives a singular perturbation problem, the solution to which 
we describe using the method of matched asymptotic expansions [55]. The solution on the domain — 00 < a; < 
comprises a two region structure, namely an outer region for —x = 0(1) and an inner region for —x = 0{ij). 

A. The outer solution 

For — a; = 0(1) we pose the regular expansions [60] 

v{x) = vo{x) + i^'^viix) + o{n^), 

Prix) = Pro{x) + iX^Prl{x) + o{tJL^) , 

Pm{x) = Pm(i{x)+p!^pml{x)+0{p!^), aS ^ ^ 0, (14) 

to give the leading order problem 

PrO+PmO = ^, -1— = , — 1 = , (15) 

K'' ax X ax x 



subject to the conditions 



At first order 0{iJ,'^) we obtain 



3x^ 

as a; ^ — oo vq ^ 1/2, pro ^ — 2~' (■'■^) 



-6axK^{pri + Pmi)(2t;o - 1)^ = - ^(1 - 2vo){2 - vo)vo^ , (17) 
dpri 4 dpml 3 

= - {VOPrl +VlPrO), — 1 = - [VoPml + WlPmo)- (18) 



ax X ax X 

with 



Eliminating p^o in (15) gives 



as X ^ — oo vi = 0(1), pri = o{l). (19) 
1. The leading order solution 



6a;2 



with pro satisfying the nonlinear first-order ODE 

dpro 24a;^/0r0 



dx {9x^ + K^Pro) 



(21) 



This last equation is scale invariant under the transformation x = ax,pro = a'^PrO for real a, suggesting the change 
of variables x = —e^jpro = x^u{y) which gives the autonomous equation 

du 2m(3 — K?u) 

with general solution 

3/4 

= (3-K2u)e«/^ 



ko 
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v„(x) 

large -x asymptotic behaviour 
small -X asymptotic befiaviour 




FIG. 1: Plot of the loading order outer solution vo{x) using (20) and (22) over a range of a; < values in the parameter case 
= feo = 1. Also shown are the large and small —x asymptotic behaviours for vo in (23) and (24). 

where fco > is an arbitrary constant of integration (and the multiphcative constants being introduced for later 
convenience). In terms of the original variables this gives the implicit solution 



2 \ 3/4 

o 2 2 o I 

dx — K PrO = o — ; — 

ko 



(22) 



with Vo then being given by (20) and PmO being determined from the first equation in (15). For reference we note the 
following limiting behaviours of the leading order outer solution, namely 



as X 



whilst 



as a; ^ 



PrO 



PmO 



2^8 \koJ 



3/4 1 / Q \ 3/2 



3 



3/4 



i-xf/' 



4^2 \ k, 



3/2 



(l_^(_,)V3), ,^3.^(_,)8/3, ^^^^1 



X — 



§(-x)«/3. 



As expected, the leading order outer solution has no dependence on a. 

We note that the outer problem has no inherent time scale, so that the scalings 

_ 7-3/2- _ - _ 7-3- _ 7,-3- 

X — Kg X, V — V, Pr — Kq Pr, Pm — ""0 Pmj 

remove the constant ko, so that without loss of generality we may consider fcg = 1- 



(23) 



(24) 



2. The first order solution 

Using the leading order solution for vo in (17) gives the simplification 

ava{3vQ + 2) 



Pml = — 



12k^x^{2vo - ly 



- Prl, 



(25) 
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and then eliminating pmi in (18) gives 



72a:4(2i;o - 1)^ 6x2 
Consequently the equation for p^i in (18) then gives 

dpri 2 , (Twg(2 - 3wo)(39wo - 2) 



av^{39vo-2) ^ 
^1 = ^o„4/o„. TT2 - Z3« Pri- (26) 



The use of an integrating factor allows us to obtain a quadrature expression for pri, however we simply note that this 
is possible rather than record it since its form will not be explicitly used. Using the behaviour (23) for vq we may 
deduce from (27) the behaviour 



as X ^ —DO, 



where fci is an arbitrary constant. The requirement (19) enforces fci = 0, thus setting the homogeneous solution terms 
to zero. Consequently from (25) and (26) we have that 



«W~-^(^fJ hx)-'^', ^^-^[f) (-^) asx^-00. (29) 
Similarly using (24) in (27) with (25) and (26), we have that 

n- 2 8ako N-4/3 2 o -2 ^(T _4 

as X ^ K pri ~ — y—x) ' , K Pml ~ — 2c7X , Vi ~ —X . (30) 



B. The inner solution 



In the outer region, the derivative of v in Eqn. (11) is neglected at leading order. This is not expected to hold 
uniformly as x 0~, since dominant balance in Eqn. (11) is given when — x = 0{p) (when using the scalings 
V = 0{l),Pr = 0{p,^^^), Pm = 0{p?) from the limiting behaviour (24)). Consequently we now consider an inner (or 
boundary layer) region at — x = 0{p), for which we introduce the scaled variables 

X = liX, V = V, Pr=II^/^Pr, pm = l^'^Pm- (31) 

The system (11)-(13) then becomes for —X = 0(1) 
dV 3(l-2y)(2-y)y 



dX 2X 

dPrr,. SVPrr,. 



dX X 

dPr WPr 



QaX{2V - If (^^{Pm + ti^'^Pr) - 3^2) , (32) 

(33) 

(34) 



dX X 
For the inner expansion, we pose 

V{X) = Vo{X) + 0{iiy% P„(X) = P„o(^) + 0(^2/3), P,(X) = P,o(^)+0(m'/'), asM^O, (35) 

to obtain the leading order equations 

dVo _ 3(l-2yo)(2-Fo)V^, 



dX 2X 

dPmO _ 3I^oPtoO 

dX ~ X ' 

dPrO ^ WpPrO 

dX X ' 



- 6aX{2Vo - i f {K^Pmo - 3X2) ^ ^gg^ 

(37) 
(38) 
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FIG. 2: Illustration of the numerical solution to the leading order inner problem (36)-(38) subject to (39) as an initial-value- 
problem solved over X £ [Xi^f , Xo] with Xi^f — — 10^,Xo = —10"'' in the parameter case = ko = 1. Both CDTT {a = 1) 
and mCDTT (ct — —1) models are considered. 



subject to the leading order outer matching conditions 

asX^-(X3 Vq^I, P„o~4(-^)', Pro - ^(-^)'/', (39) 

which follow from (24). It is noteworthy that Pro does not enter (36), illustrating that it does not affect the leading 
order behaviour of V{X) in this inner region, only its correction terms. Again ko — 1 may be considered without loss 
of generality since the scalings 

X = X, Vq = Vq, Pro ~ koPj-o, Pm — -PmOi 

remove the constant ko from the leading order inner problem. 

An eigenmode analysis given in Appendix A about the far-field asymptotic behaviour (39) for the system (36)- 
(38), confirms that (39) impose the necessary three conditions on the nonlinear third order system. Consequently 
this leading order inner problem may be solved as an initial- value-problem (IVP). Using (36)-(38), a more accurate 
far-field expansion consistent with (39) is 

asX->-oo Vo^- + —X-\ p„^^^ —{-Xfe-^, Pro ^{~Xf/^e-^ . (40) 

Using these asymptotic expansions, Eqn (36) becomes 

dVa 40 a 



dX 9 X 



5 



(41) 



and we see that the slope of the asymptotic behaviour depends on a. Numerical solution of (36)-(38) subject to (40) 
was obtained using MATLAB's IVP solver odelSs, using solver tolerances of AbsTol=RelTol= 10~^^. The system 
was solved over the interval X S [Xinf,Xo] where Xinf = — 10^, ATq — —10"'^, the condition (40) being imposed at 
X = Xinf- The parameter values = kg = 1 where used. Figure 2 illustrates the behaviour of Vq, for the cases 
cr = ±1. 

A local analysis near Vq — 2 for the system (36)-(38) gives the behaviours 

486cr A , c 729crX^ , Q 486tTX^ 

asX^O- Vo^2-^X\ P„o-AmO^'e P.o ~ Ao^'e —, (42) 

5 

for constants v4„io,^rO (which can be determined numerically A„ia ~ 98,^^0 ~ 105 in the case = fco = 1 and 
a = l). 

For —X <C 1, (36) simplifies to 

dVo _ 3(l~2Uo)(2-Uo)yo 

dX ~ 2X ' ^ ' 
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■Mjyj _ 3(14il<;)(2-\<S)V 
dX 2X 
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FIG. 3: The phase plane of possible solution branches of the leading order inner problem, Eqn. (43), which is independent of 
o". We have included negative values of v for completeness, although we never consider such values cosmologically. Notice that 
the solutions are symmetric around X —X. 

a phase plane of possible solution branches of which being given in Figure 3. The numerical solution of Figure 2 enters 
the stable equilibrium point Vq ~2 along one of the branches shown in {— oo < X < 0, 1/2 < Vq < 2}. 

The dependence of (41) on a is crucial for the inner solution. For the CDTT model, the positive slope ensures that 
the solution moves from Vb ~ 2/3 {X — > — oo) to the final solution Vb ~ 2 (X ^ 0~). However, when a — —1, the 
asymptotic behaviour around Vb ~ 2/3 has negative slope and the solution rapidly approaches the solution, Vb = 1/2. 
We believe this is the "recent radiation epoch" behaviour found in [40, 50], which is not cosmologically viable. As 
X ^ 0~, the numerical solution then follows that of Figure 3 and Vb 2. In fact, the a = —1 (CDTT) model cannot 
accomodate the outer solution of Vq ~ 2/3, but rather must match = 1/2 for the leading order outer solution. As 
such, then, we do not obtain a viable asymptotic solution in the a = —1 case. In the remaining sections of this paper, 
we only consider the mCDTT model, which has a cosmologically viable solution. 

C. Matching and composite expansion of the mCDTT model 

For completeness we now illustrate matching between the two term outer solution and the one term inner solution 
for the mCDTT model, before considering a composite asymptotic expansion valid for all a; < 0. 

Taking two terms of the outer expansion (14) and using the behaviours (24) and (30), we have for small —x that 

v{x) = v,{x) + ^,^v,{x) + o(/) §^(-2;)2/3 + 

Pr{x) = pM + p^^p,,{x) + o{p^) ^ %{^xf/^ //^(_^)-4/3^ 

Pm{x) = Pmo{x) + fJ.'^Pmlix) + o{p^) ^ ^X^ - ^{-xf^^ - p^^X^"^, 

which in terms of inner variables (31) becomes 

P.(X)^^(-X)«/3-|^(-X)-4/3, 

P„AX) - ^X' - + 0(^2/3)^ (44) 
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Composite expansion 



comp 




-10^" -10^ -10° -10"^ -10"'° 



X 



FIG. 4: Illustration of the leading order composite expansion when a — 1 (mCDTT model) 
for V given in (46) for selected fi values. 



Taking the inner expansion (35) and expanding (40) for large gives 

V{X) = Vo{X) + 0(p'^') = 1 + ^X-^ + o{X-^) + 0(/.2/3), 

UX) = Pr,{x) + o(^^/3) = ^i-xf^' (i - gji) + o{{-x)-'/') + 0{^.y% 

Pm{X) = PMX) + 0{i/'') = ^X' [l ^) + o{X-') + 0{^,^'% (45) 

these expressions agreeing with those in (44) to the order of the terms retained (i.e. the first two terms). 
A composite expansion using the leading order terms of the outer and inner expansions is given by 

Vcompix) = vo(a;) + Vo - ^, (46) 

Pr-comp(a;) = pM + Pro ( - ) " ^{-xf'"^, (47) 
Pmcomp(a:) = Pmo{x) + PmO ( " ) ^ (48) 

where the last terms in each expression are the leading order terms in the overlap region between —x = 0(1) and 
—X = 0{fi). Such expressions are now uniformly valid throughout the interval — oo < a; < 0. A plot of Wcomp(a;) is 
given in Figure 4 for three selected values = 10^^, 10"'^, 10^^ illustrating the emergence of a two step profile as n 
decreases. 



V. COSMOLOGICAL IMPLICATIONS 



Since there is no observation that can probe the future expansion of the universe, we are restricted to testing the 
validity of the outer solution, which describes the past: 

v{x) = vo{x)+fi'^vi{x)+o{fi'^). (49) 

This can be done in several ways, but we restrict ourselves to just two: the requirement that we have a long enough 
period of matter domination to satisfy large scale structure and that we satisfy the constraints of Big Bang Nu- 
cleosynthesis (BBN). We restrict ourselves to the mCDTT model, which we have seen leads to a standard cosmic 
history. 
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We define an "initial" epoch, given by xq, at which the ratio of radiation to matter is set: 

ro = (50) 

We can relate the parameter fco to this ratio, 

fco = 3(ro + l)*/'(-a;o)-5 (51) 

Setting ko = = l (as in Section IV), we find matter-radiation equality (rg = 1) at xq = — 12-\/3 « —20.78, as can 
be seen in Figure 4. 

In comparison to data, we must convert the "timescale" from H (x) to redshift, z: 

where ao is the scale factor at present times. We define this to be the transition time from v = | (such that the 
universe is just accelerating today) and acts as a normaUsation. We note that 

V 

dlna = Hdt = dx. (53) 



A. Sufficient Matter Domination 



The redshift of equality, Zgq, is found from the first peak in the Cosmic Microwave Background (CMB) anisotropy 
spectrum. From the WMAP first year data [57], Zcq = 3454^332 (WMAP data only). In order to ensure a long enough 
epoch of matter domination, we require that v{x) « 2/3 up to this redshift. 

We may relate the redshift of matter domination back to the free parameter ko. We may use Eqn. (53) and require 
that, at the present time, x{0) = —Hq = — /x and Inao = 0. Hq is the present value of the Hubble parameter. 
Assuming that t; = |, we find the relation: 

a: = a:o{^f" = -^{l + zf\ (54) 

At radiation-domination equality, Xeq ~ — 2 x 10^/i. Using Eqn. (51) and defining fco at Xeq, we find that, to ensure a 
long enough period of matter domination: 

fco w 0.002/i-2/3. (55) 

Despite appearances, we note that this solution agrees with the analysis of [52], although we quantify the period 
of matter domination and find a sufficiently long era. We note that, during this matter epoch, due to the asymptotic 
solution, flm 1, but never reaches (nor sits at) unity. In Amendola et al's notation, the solution undergoes a 
transition from the critical point, P5 to Pq, where — > 0~, me = —2 and mg g > —1. This transition is allowed 

according to [52]. At the transition, such variables (m etc.) become infinite, due to the factor of 1/(1 + /). 



B. Big Bang Nucleosynthesis 

One of the most constraining tests of non-standard cosmological expansion rates comes from the earliest epoch we 
can test: Big Bang Nucleosynthesis. In an early radiation-dominated epoch, the expansion rate is determined by the 
energy density of the relativistic particles. Due to this, tight experimental constraints exist. (For a review, see [58]). 
Firstly, we note that at the time of BBN, ^;bbn = 5.9 x 10^, and using Eqn. (54), 

^=-^ = (5.9x10^)'/'. (56) 
One parameterisation of the BBN constraint is given as 

Pr ^ p'r = Pr + PX = S"^ Pr 
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Composite expansion 




z+1 



FIG. 5: The composite solution, Wcomp shown against redshift, z, in the mCDTT model (cr = 1) 
when /X = 10~^Mpi. In order to obtain sufficient matter domination, we require fco = ^."iM^^^^ . 



where px denotes the extra component, in our context due to gravity. The constraint parameter, S, can be written 
in terms of the addition of an extra neutrino species, 



5^ = 1 + 7/43A7V^ < 1.16, 



(57) 



where, in the last term, we have assumed that AiV,^ < 1 [59]. In the mCDTT model, with the asymptotic solution 
found in Section IV, we find 



S-^-l = 



PrO 



8fco 

'9k^ 



•J-^BBN 



n -1 



—fi ko (-xbbn) 



Substituting the value of kg required for sufhcient matter domination, Eqn. (55), we find 

5^ - 1 w -5.9 X 10""* 
which is obviously well within present constraints. 



10/3 

— — ) «-8.25xl0~^^ 

-XBBN 



VI. GENERALISATION 

It must be acknowledged that the mCDTT model may ruled out by observation within the solar system and 
therefore is perhaps not the best model to consider. We note, however, that it provides a useful example of the 
technique and illustrates at least one cosmologically viable f{R) model. We now wish to consider other models for 
which this technique is valid, so to identify a class of models which might work fully. 

We firstly note that for any model with R + f{R), the system can be termed as in Eqn. (9) and therefore the leading 
order solutions, given by Eqns. (22-24), are always true. This is to be expected, since we want to pick up the Einstein 
solution when fi 0. Due to this solution, we also expect to have the same inner scaling for v, pr and p™, as given 
in Eqn. (31): 



V, Pr=e^/^Pr 



1 rr 
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where we have substituted e = /i^ to account for a more general parameterisation. With these scalings, the inner 
problem for a general f{R) theory can be written: 

dV _ V^jV-l) fR / 2^(2^-1) [3x''-K^P^ + fiy^Pr)]V^ 

dx 6x3 ^1/2 36x^efnR ^ 18x^ e^/^fRR 

dPm ^ SVPm 

dX X ' 

dPr WPr 



dX X 

The relevant scaling is suggested from e^/^/ = 0{1) (although this is not the only condition). For inverse power-law 
models, 

~f{R) = AR-^ (59) 
with positive n, we find the scaling x = e^/^"X so that 

dv_ _ y(y-i)(2y-i) v{2v-if 2^(2^-1) 2(6")y^-"(2y - ly+^x^ri-i [3^2 _ ^2 ^ /jV^p^^j 

dX ~ {n + l)X n(n+l)X ^ X n(n + l)yl 

y(2y- l)(n + 2) 



n(n + l)X 



^_ (n+l)(2n+l) 



(n + 2) 



2(6")V-l-"(2y - 1)"+2X2"-1 [3X^ - /t^ {Pm + l^'^'^Pr)] 

n{n + l)A 



We thereby note that as X — > — cxd, the leading order term is the last term (proportional to X^"~^) and we expect to 
recover the Einstein solution. As X ^ 0^, the solution is dominated by the first terms. We expect dominant balance 
when — a; = 0(e^/^") and therefore sufficient matter domination requires 

fco = 0.002e-^/''". (60) 
The (accelerating) asymptotic solution as X — > 0~, which we denote by 14tt, can be calculated from 

_ (n+l)(2n+l) 

• ^^^^ 

For n = 1, we find T4tt = 2 as seen in Figure 4. We strongly note, however, that the sign of the coefficient. A, seems 
crucial if one wishes to avoid the apparent ^MDE solution of [50] . In each model, a full matching analysis is required. 



VII. CONCLUSION 



A modification of gravity has been suggested in the form of f{R) theories. However, until now, it has been believed 
that these cannot produce a standard era of matter domination. By studying the asymptotic behaviour of a class of 
f{R) models, we show that a sufficiently-long period of matter domination can be achieved. 

We consider one specific example, the (m)CDTT model, R + aji^ /R with a = ±1, and explicitly calculate the 
asymptotic behaviour. In terms of the CDDETT variables, where x = —H, we find that the system gives a singular 
perturbation problem, comprising a two region structure: an outer region —x = 0(1) and an inner region x = 0{ijl). 
We show that the mCDTT solution recovers Einstein cosmology when —x = 0(1) (i.e. in the past when H ^ fi) 
and matches an accelerated cosmology when —x = 0{fi) (i.e. today when H ^ (x ^ Hq). We directly relate this 
solution to cosmological variables and find the condition for which we obtain sufficient matter domination to satisfy 
large scale structure data. We also show that the perturbation causes no significant effect at the time of Big Bang 
Nucleosynthesis . 

We generalise our results to include the class of inverse power-law functions. Although these are also ruled out by 
solar system tests, we aim to exhibit the generality of our previous results, which could be valid for other functions. 
We find that the original CDTT model cannot match to the outer ?; ~ 2/3 (matter dominated) solution and the 
"recent radiation" behaviour of [40, 50] is found. Generally, we find that for models with R ± eAR^^, standard 
cosmological evolution is expected until —x = 0(e^/^"), which can be coincident with the present time. After this 
point, an accelerated solution is found. This solution depends strongly on the sign of A. 

Our main result is that, for a model with lagrangian R + .f{R), we expect the leading order outer (past-time) 
solution to be given by Eqn. 22, the Einstein solution. The next order terms will be proportional to the perturbative 
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parameter in f{R). The leading order inner (future) solution can be found using the system in Eqn. (58). We therefore 
expect the solutions to match across the boundary and the universe should find the accelerating attractor, given by 
Eqn. (61). Whether the "recent radiation epoch" is found (instead of a matter epoch) seems to depend on the sign of 
the corrective term. 

It is the authors' belief that this perturbation method can be applied to even more general functions of higher 
curvature invariants, for example f{R, P, Q), where P = R^^R^'^ and Q = RapjsR"^''^ , such as the CDDETT model. 
This is the focus of future work. With solutions such as we have shown here, observational constraints can be applied 
to a much larger class of modified gravity models than is possible at present. 



APPENDIX A: EIGENMODE ANALYSIS (MCDTT MODEL) 



The number of degrees of freedom contained within the asymptotic behaviour (39) may be determined through an 
eigenmode analysis. We consideration a perturbation about the far-field behaviour in the form 



Vo{X) 



■dv{x), PMX) 



+ SP^iX), Pro{X) ~ ■^(-X)8/3 + SPriX), 



Keeping terms at 0{S), (36)-(38) give the linearised equations 



X 



dX 

dPr 

Ix 



as X -oo.(Al) 

(A2) 
(A3) 
(A4) 



Eliminating P^ from (A2) and (A3) gives the transformed Bessel equation 



X' 



,fv_ 

dX2 



+ 4X(1 - + 2{X^ - 6)t> = 



with solutions 



r {-xrMxyV2), 

^ {-XrY4Xy^/2), 



iV57- 24(7 + 16cr2 



where 

m — 2(7 — 3/2, ly 

where J,y{z),Yi,{z) are standard first kind Bessel functions of order v. Their large —X behaviours are 

(_X)'cos(fJ-57r), 



V 



where 



(-X)'sin(^-5.) 



; = 2a-5/2, S 



as X —00, 



(21/ + 1) 



(to within multiplicative constants) and consequently the three linearly independent eigenmodes are 



V ~ {-Xy COS - Sir 

2^(_X)'+V3 sin - Sn) 



y~-(-X)'sin(^-5.) 

_|_(_X)'cos(^-5.) 



Pr.. 

p 2\/2fci 

It 



^(-X)Wcos(^-57r) 



> ; 



(A5) 



(A6) 



V = 

Pm=0 

P^ = (-X)8/3 



.(A7) 



The third mode in (A7) corresponds to small changes in ko and is the only mode consistent with the far-field behaviour 

(39) . This implies that there is only one degree of freedom in the asymptotic behaviour (39) which is associated with 
ko- Thus if ko is specified then the far-field behaviour imposes three conditions on the system (36)-(38). 
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